Silencing of GhORP_A02 enhances drought tolerance in Gossypium hirsutum

Background ORP (Oxysterol-binding protein-related proteins) genes play a role in lipid metabolism, vesicular transferring and signaling, and non-vesicular sterol transport. However, no systematic identification and analysis of ORP genes have been reported in cotton. Result In this study, we identified 14, 14, 7, and 7 ORP genes in G. hirsutum, G. barbadense, G. arboreum, and G. raimondii, respectively. Phylogenetic analysis showed that all ORP genes could be classified into four groups. Gene structure and conserved motif analysis suggest that the function of this gene family was conserved. The Ka/Ks analysis showed that this gene family was exposed to purifying selection during evolution. Transcriptome data showed that four ORP genes, especially GhORP_A02, were induced by abiotic stress treatment. The cis-acting elements in the ORP promoters were responsive to phytohormones and various abiotic stresses. The silenced plants of GhORP_A02 were more sensitive to drought stress when compared to control. Conclusion The major finding of this study shed light on the potential role of ORP genes in abiotic stress and provided a fundamental resource for further analysis in cotton. Supplementary Information The online version contains supplementary material available at 10.1186/s12864-022-09099-y.

ORP (oxysterol-binding protein-related proteins) genes play a key role in lipid metabolism, vesicular transferring and signaling, and non-vesicular sterol transport [41]. Previous studies about ORP genes in Arabidopsis, soybean, and petunia have also demonstrated its significant role in biotic stress, abiotic stress [31,43,45]. The Arabidopsis genome encodes 12 ORP genes, and the rice genome encodes six ORP genes [48]. Although the ORP gene in plants has been cloned, there are few studies focused on their functions. In Arabidopsis, ORP3a, located the endoplasmic reticulum, interacts with VAP33 family member PVA12 [43]. In Petunia inflata, PiORP1 participates in pollen growth and development by interacting with PRK1 receptor kinase on the plasma membrane of a hybrid pollen tube [45]. In soybean, the expression of GmOSBP was inhibited by salt stress but induced in aging leaves, indicating that GmOSBP may be involved in stress response and the cell aging process [31]. In this study, we performed genome-wide identification and investigated phylogenetic relationships, gene structure, conserved domains, gene duplication events and expression files of ORP genes. Our study may be useful for the future molecular and biological function of the ORP gene family in cotton.

ORP gene identification in cotton species
The conserved domain PF15413 of ORP genes was obtained from PFAM (http:// pfam. xfam. org) and used as a query sequence to retrieve the ORP genes in four cotton species by Hmmer 3.0 (http:// hmmer. org), and the identity of the ORPs genes was analyzed by SMART (http:// smart. embl-heide lberg. de/ smart/). The physical and chemical characteristics of ORP proteins, including molecular weight, protein length, molecular charge, isoelectric point, and grand average of hydropathy, were obtained from CottonFGD (https:// cotto nfgd. net/).

Chromosomal mapping
We used the GFF3 files of the ORPs genes downloaded from CottonFGD to find the distribution on all chromosomes. TBtools software (version 1.098685) was then used to visualize the gene's location on chromosomes [8].

Phylogenetic tree and collinearity analysis
The full-length protein sequences of ORP genes from Gossypium were downloaded and aligned using ClustalW with default settings. The phylogenetic tree was constructed using the neighbor-joining method in MEGA 6 with default parameters and 1000-bootstrap replicates (http:// www. megas oftwa re. net/). The protein sequences of Gossypium hirsutum (G. hirsutum) have been searched along the protein databases of Gossypium arboreum (G.

arboreum), Gossypium barbasense (G. barbasense) and
Gossypium raimondii (G. raimondii) by BlastP to identify homologous genes and hits with E-values of 1.0E -5 and similarity of 90% were considered noteworthy. TBtools program was used to create the collinearity analysis using the GFF3 file, linked file, and gene IDs. Collinearity analysis was performed among three cotton species (G. hirsutum, G. arboreum, and G. raimondii) using Circle gene viewer in TBtools software to determine collinear gene pairs. Coding and protein sequences of all homolog genes were used to calculate the Ka/Ks (Non-synonymous substitution-rate/Synonymous substitution rate) value by TBtools [47].

Expression profile analysis of ORP gene family
FPKM values (fragments per kilobase of exon per million fragments mapped) of ORP genes were downloaded from CottonFGD. We analyzed the expression profiles of ORP genes under different stress treatments, which included PEG, salt, heat and cold treatments.

Virus-induced gene silencing (VIGS)
For virus-induced gene silencing, the cotton variety H117 was employed. H177 was developed by the Institute of Cotton Research Anyang of the Chinese Academy of Agricultural Sciences. This Variety was chosen because it is particularly susceptible to many environmental stresses, including drought. A 306-bp fragment of GH_A02G0809 was amplified from G. hirsutum acc. TM-1 with gene-specific primers. The PCR product was then digested with Spe I and Acs I and cloned into Spe I-Acs I -Cut pCLCrVA. The fusion vector was named pCLCrVA: GhORP_A02 and transformed into Agrobacterium tumefaciens strain LBA4404. The control vector pCLCrVA, pCLCrVA: GhORP_A02 and positive vector pCLCrVA: PDS were mixed with pCLCrVB at a 1:1 ratio [19]. The mixed Agrobacterium tumefaciens solutions were injected into the ten-day-old cotton cotyledons on the abaxial side with a needle-free syringe. The plants were placed at room temperature in the dark overnight and grew at 23 °C with a 16 h / 8 h light/dark cycle. Agrobacterium infection was carried out three times with 30 plants for each vector. The primers for VIGS vector construction are listed in Table S1. Wild type and the plants injected with pCLCrVA empty control and pCLCrVA: GhORP_A02 were subjected to drought treatment after four weeks. Drought treatments of the seedlings were irrigated with 15% PEG6000, while control plants were irrigated with 1/2 MS nutrient solution.

RNA extraction and quantitative real-time PCR (qRT-PCR) analysis
Total RNA was extracted from fresh leaves and roots using TRIzol ® Plus RNA Purification Kit (Invitrogen, CA) based on the manufacturer's instructions. Approximately 1 µg RNA was reversely synthesized into cDNA using the iScriptTM Synthesis Kit (Quanta BioSciences, MD). The qRT-PCR was carried out in an Eppendorf real-time PCR equipment using a 5 µl cDNA template (diluted 1/100), 5 µl primers (2.4 M), and 10 µl SYBR green mixture (Promega, Madison, WI). Histone 3 was used as the internal control, and the relative expression levels of the ORP gene were calculated by the 2 −ΔΔCt method [35].

Physiological analysis
Physiological parameters, including ion leakage, chlorophyll content, excised leaf water loss, and relative leaf water content, were determined after 10 days of drought treatment. Wild type, silenced, and control plants (ten plants for each) were harvested after drought stress for oxidant and antioxidant concentration analysis. The H 2 O 2 content, peroxidase (POD), malondialdehyde (MDA) and catalase (CAT) were determined by using the corresponding ROS content reagent kits and enzyme activity kit (Solarbio, China) according to the manufacturer's instructions. The experiment was repeated three times.

Genome-wide identification and chromosomal locations of the cotton ORP genes
To identify all ORP genes in two allotetraploid cotton, G. hirsutum (AD 1 ), G. barbadense (AD 2 ) and its two diploid ancestors G. arboreum (AA) and G. raimondii (DD), we used conserved domain Pfam 15,413 to retrieve ORP genes and identified 42 ORP genes. G. hirsutum, G. barbadense, G. arboreum, and G. raimondii have 14, 14, seven, and seven ORP genes, respectively. The ORP gene ID and predicted protein properties and subcellular locations are listed in Table 1. Variable distribution of ORP genes on chromosomes across all four cotton species was observed (Fig. 1). In G. hirsutum and G. barbadense, ORP genes were uniformly distributed on the At and Dt chromosome. In G. hirsutum and G. barbadense, 14 ORP genes were located on chromosomes A02, A03, A05, A06, A09, D02, D03, D05, D06 and D09. Two ORP genes were located on chromosomes A03, A05, D02 and D05. In G. arboreum, seven ORP genes were located on chromosomes A01, A03, A05, A06 and A09. In G. raimondii, seven ORP genes were located on chromosomes D03, D05, D06, D09 and D10.

Phylogenetic and collinearity analysis of cotton ORP genes
According to the phylogenetic analyses, all the ORP genes could be classified into four clades ( Fig. 2A). Four, two, two, and six ORP genes from G. hirsutum were classified into Group I to Group IV. Both in G. arboreum and G. raimondii, Group I, II, III and IV have two, one, one, and three ORP genes. To analyze the evolution of the ORP genes from diploid to tetraploid species, collinearity analysis was performed among three cotton species (G. hirsutum, G. arboreum, and G. raimondii). There were seven, seven, and seven orthologous gene pairs between the A and D genomes, the At subgenome and the A genome, and the Dt subgenome and the D genome (Fig. 2B). The number and relatedness of ORP genes in the three species suggested that ORP genes were not lost during G. hirsutum speciation.

Evolution of ORP genes in Gossypium species
Natural selection has no effect on gene's Ka/Ks values during the evolutionary trend, but Ka/Ks > 1, Ks/ Ka = 1, or Ka/Ks < 1, the Ka/Ks value indicates positive, neutral, or negative selection, respectively [55]. Similar results were found in the distributions of Ka, Ks, and Ka/Ks among homologous pairs of Gossypium species. The Ka/Ks ratio of most orthologous gene pairs was less than one, indicating purifying selection during evolution resulting in limiting the functional divergence after duplications and polyploidization of ORP genes (Table  S2). Only two orthologous gene pairs (GH_A02G0809 and Ga3G0877, GH_D02G0824 and Gorai.005G091400) have Ka/Ks ratio exceeding one, which implies these gene pairs underwent positive selection and had relatively rapid evolution rate.

Gene structure and motif identification of ORP proteins
The gene structure of ORP genes was analyzed according to the annotation files (Fig. 3A-D). Most ORP genes have 8-10 exons and only six genes in four cotton species have two or three exons. Genes classified into the same evolutionary branch have conserved gene structure patterns in terms of exon number and exon length. The MEME search identified ten conserved motifs in ORP genes, ranging from 300 to 2100 amino acids (Fig. 3E-H). The conserved motif numbers in different genes varied from 3 to 10. Motifs 2, 5, and 7 were conserved in all ORP genes in G. hirsutum. In G. barbadense, the conserved motifs were 4, 5, and 7, while in G. arboreum and G. raimondii, the conserved motifs were 4, 6, and 7.

Identification and analysis of cis-acting elements
Cis-acting regulatory elements play a key role in molecular switches that control a dynamic gene activity network that initiates many biological processes, such as hormone responses, developmental processes,   G. hirsutum, B G. raimondii, C G. arboreum, D G. barbadense. Based on their genome, each species' chromosomal location was plotted and abiotic stress responses [37]. MBS (drought inducibility), ABRE (abscisic acid-responsive), and TC-rich repeats (repeat cis-actin), which are involved in defense stress response and drought stress, could be found in all four cotton species ( Figure S1). CAT-Box (meristematic cell expression), GARE-motif (Gibberellin responsive), and TGA-elements (auxin-responsive element), which are involved in the germination and regeneration stage, could also be found in all four cotton species.

Expression profiles of ORP genes in Gossypium hirsutum
The raw RNA-seq data of the 14 ORP genes in G. hirsutum were normalized to log 2 (FPKM), and the heatmap of the expression is presented as Fig. 4. Four genes, including GH_A02G0809, GH_D02G0824, GH_A09G2118, and GH_D05G2343, were up-regulated by the abiotic stress treatment. Especially, GH_A02G0809 (GhORP_A02) expression was induced significantly by PEG treatment. We further performed experiments to characterize the function of GhORP_A02 in drought stress.

Virus-induced gene silencing of GhORP_A02 in cotton show significant sensitivity to drought
The method for gene silencing through virus-induced was used to analyze the role of the GhORP_A02 in drought tolerance. Gossypium hirsutum acc. H177 was infected with three vectors, including pCLCrVA: PDS (positive control), pCLCrVA (negative control), and pCLCrVA: GhORP_A02. Ten days after infection, the indicator pCLCrVA: PDS showed albino color, the control plant showed a normal color without visible change, and the pCLCrVA: GhORP_A02 plants showed complete shrinkage of the leaves, which indicates that VIGS was successful (Fig. 5A). qRT-PCR was used to analyze the expression level of GhORP_A02 in silenced plants, and the result showed that the infected plant (pCLCrVA: GhORP_A02) showed a lower expression level than the control plant (Fig. 5B). The physiological analysis includes ion leakage, chlorophyll contents, excised leaf water loss, and relative leaf water content was done in silenced and controlled plants with and without drought treatment. The relative ion leakage level of the silenced plants increased by 20% compared to the control. The chlorophyll contents of the silenced plant were significantly lower in comparison to control plants. While in excised leaf water loss, the silenced plants lost more water than the control plant. The relative leaf water content of the silenced plant didn't show a significant difference when compared to the control under drought conditions (Fig. 6A-D). Determination of antioxidant (CAT and POD) and oxidant (MDA and H 2 O 2 ) enzyme concentration levels were analyzed in both control and silenced plants under drought conditions. There is a significant increase in the concentration of antioxidants and a decrease in the concentration of H 2 O 2 in silenced plants compared with their respective control (Fig. 6E-H).

Discussion
Drought is one of the most significant abiotic stresses, resulting in considerable yield losses in cotton [23,40]. Plants have evolved self-defense systems to deal with abiotic stresses, which involves the transcription of stressrelated genes [39]. Genetic enhancement of drought tolerance hinges on identifying genes related to drought tolerance [53]. In earlier research, drought-responsive genes were identified in many species like rice, peanut, soybean, wheat, maize and cotton [9,21,36,44,54]. Oxysterol-binding protein (ORP) and its homologs constitute a protein family in many eukaryotes, from yeast to humans, which are involved in cellular lipid metabolism, vesicle transport and signal transduction [51]. Recent studies have demonstrated that the ORP gene family was stress-responsive in various plants [31,48,43,45]. Our current research has demonstrated the function of GhORP_A02 in drought stress response. This study used the protein domain PFAM 15,413 to retrieve ORP genes in the four cotton species, and G. hirsutum, G. barbadense, G. raimondii, and G. arboreum encoding 14, 14, seven, and seven ORP genes, respectively. In previous research, 12 and six ORP genes were identified in Arabidopsis and rice, respectively [48]. Gene structure and phylogenetic tree analyses indicated that all GhORP genes, classified into one group, have a similar gene structure. The evolution analysis of the ORP gene in four cotton species shows no negative selection across all the species.
Conserved domains correspond conformational changes due to binding [28,50]. Domain rearrangement and recombination, which typically occurs due to gene duplication and fission or fusion events, are used to develop new protein functions [38]. In this study, we identified two conserved motifs, motifs 5 and 7 in G. hirsutum and G. barbadense, while motif 7 was conserved in all four Gossypium species. Subcellular localization and the transcription of a gene under stress are powerful mechanisms to explain its biological function. ORP genes have been identified in Arabidopsis, soybean, rice, and Petunia and found to be located in plasm membrane, nucleus, and endoplasmic reticulum [31,48,43]. According to the subcellular localization prediction result, most ORP genes are located in the nucleus. Expression analysis of GhORP genes under different stress showed four genes were strongly induced by cold and drought stress from 1 to 24 h. A similar result was reported in soybean and Arabidopsis. For example, GmOSBP was induced by salt stress, and AtORP4A and AtORP4B were induced by drought stress [31,48]. Expression analysis showed that four GhORP genes, especially GhORP_A02, were significantly up-regulated after drought stress, and we consider this gene as the candidate gene for drought stress response.
To investigate the function of GhORP_A02, we silenced this gene by VIGS. Resultantly, silenced plants were more sensitive to PEG treatment than control. Environmental challenges such as drought, salt and temperature cause a redox imbalance in plant cells, which rises the total rate of metabolism and finally upregulates H 2 O 2 production [17]. There is still no relative study about the mechanism of ORP proteins to cope with abiotic stress. How ORP proteins are involved in stress response is largely unknown. When plants are subjected to abiotic stress, membrane proteins degrade, and comparative conductivity and MDA are significantly elevated [12,25]. In this study, electrolytes in silenced plants (pCLCrVA: GhORP_A02) increase significantly under drought stress compared to control plants. Both chlorophyll content and relative water contents decrease significantly in silenced plants, and this is in agreement with many previous research findings, which indicated that plants tend to close stomata to avoid water loss and decrease photosynthesis in drought conditions [22,29,46]. The reactive oxygen system produces substances such as POD and CAT, which are accompanied by an increase in reactive oxygen to limit and regulate the damage of reactive oxygen to plants, but also as a signal molecule to activate the plant body to respond to the external adverse environment [2,3]. Our present work showed that Both CAT and POD decreased in the pCLCrVA: GhORP_A02 plants when compared to control plants, and this signifies the signaling role of the ORP gene in enzymatic activity in cotton, which is consistent with previous research [33,42].